.libPaths()
# delete data in global environment:
rm(list = ls())
# load standard packages:
library("stringr", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("reshape", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("ez", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("ggplot2", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("plyr", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("nortest", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("car", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("ltm", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("lsr", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("afex", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("nlme", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("ggthemes", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("robustHD", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("cowplot", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("tidyr", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")
library("lattice", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")                                
library("plyr", lib.loc="/Library/Frameworks/R.framework/Versions/3.2/Resources/library")

# load some necessary functions, for within-subject, repeated measure designs
load("/Users/XXX/effect.size.Rdata")
load("/Users/XXX/mean.sd.Rdata")
load("/Users/XXX/min.mean.sd.max.Rdata")
load("/Users/XXX/summarySE.Rdata")
load("/Users/XXX/summarySEwithin.Rdata")
load("/Users/XXX/normDataWithin.RData")
load("/Users/XXX/rcontrast.RData")

#read csv file:
getwd()
setwd("XXX")
data_for_r = read.csv("data_time_exp_future_banana.csv")

#delete NAs
completeFun <- function(data, desiredCols) {
  completeVec <- complete.cases(data[, desiredCols])
  return(data[completeVec, ])
}
data_for_r<-completeFun(data_for_r, "rt")

# only participants that have more than 50% correct answers per condition & response:
data_for_r<-subset(data_for_r, 
                   participant==1 |
                     participant==2 |
                     participant==5 |
                     participant==6 |
                     participant==7 |
                     participant==8 |
                     participant==9 |
                     participant==10 |
                     participant==15 |
                     participant==16 |
                     participant==17 |
                     participant==18 |
                     participant==19 |
                     participant==21 |
                     participant==22 |
                     participant==23 |
                     participant==24 |
                     participant==26 |
                     participant==27 |
                     participant==31 |
                     participant==32 |
                     participant==33 |
                     participant== 34|
                     participant== 36|
                     participant== 39|
                     participant== 40|
                     participant== 45|
                     participant== 47|
                     participant== 48|
                     participant== 49|
                     participant== 52|
                     participant== 53|
                     participant== 54|
                     participant== 55|
                     participant== 56|
                     participant== 57)

# response time: subtract 666ms (= stimulus length)
data_for_r$rt_offset<-data_for_r$rt-666 # -> total amount of items: 4320
# only Responses that were given before stimulus onset:
data_for_r<-subset(data_for_r, rt_offset>0) # 12 deleted
data_for_r<-subset(data_for_r, rt_offset<6000) # 26 deleted
# 38 (out of 4320) outliers = .8 %

#select only 'test blocks' (without learning blocks) per condition
data_for_r$test<-grepl(pattern = "(test)", data_for_r$blockcode) 
data_for_r<-subset(data_for_r, test=="TRUE")  

# rename conditions:
data_for_r<-mutate(data_for_r,blockcode=as.character(blockcode))
data_for_r<-mutate(data_for_r,blockcode=sapply(strsplit(data_for_r$blockcode, split='_', fixed=TRUE),function(x) (x[2])))
data_for_r$blockcode=gsub(data_for_r$blockcode, pattern = "vorwaerts", replacement = "forward")
data_for_r$blockcode=gsub(data_for_r$blockcode, pattern = "rueckwaerts", replacement = "backward")
data_for_r$blockcode=gsub(data_for_r$blockcode, pattern = "stehen", replacement = "standing")
data_for_r$condition <-data_for_r$blockcode

# only trials were an answer (correct or incorrect) was given
data_correct_and_incorrect<-subset(data_for_r, correct=="correct" | correct=="incorrect")  

# only correct answers:
data_correct<-subset(data_for_r, correct=="correct") 

#select relevant columns:
data_correct<-data_correct[c("participant", "Response", "condition", "correct", "rt_offset", "order","trialnum", "blocknum", "age")]

# count number of correct and incorrect items:
data_correct_count<-count(data_correct, c("correct", "participant", "condition", "Response", "order"))

# descriptiv (summary for within-participant designs)
# per response:
dataR <- summarySEwithin(data_correct, measurevar="rt_offset", withinvars=c("Response"), idvar="participant")
dataR # -> rt of "past" items was around 50ms faster than "future" items
# per condition:
dataCon <- summarySEwithin(data_correct, measurevar="rt_offset", withinvars=c("condition"), idvar="participant")
dataCon # -> mean rt of in the backward condition was around 50ms slower than in forward condition and around 100ms slower than in the standing condition
# per response & condition
datac <- summarySEwithin(data_correct, measurevar="rt_offset", withinvars=c("condition", "Response"), idvar="participant")
datac
# -> "future" items during backward walking were around 150ms slower compared to "past" items. 
# SDs from 359ms (standing - future) to 452ms (backward - future)

## PLOTS ##
# with outliers # 
ggplot(datac, aes(x=condition, y=rt_offset, fill=Response)) +
  geom_bar(position=position_dodge(.9), colour="black", stat="identity") +
  geom_errorbar(position=position_dodge(.9), width=.25, aes(ymin=rt_offset-ci, ymax=rt_offset+ci)) +
  theme_classic() +
  theme(legend.position="top") +
  xlab("\nCondition") +
  ylab("Response Time (ms)\n") +
  labs("d") +
  coord_cartesian(ylim=c(2500,1000)) +
  scale_fill_manual(values=c("#CCCCCC","#FFFFFF")) +
  theme_bw() +
  theme(text = element_text(size=30))+
   geom_hline(yintercept=38) 


## STATISTIK ##

# Setting contrasts
# We compare each condition against the control condition
data_correct$condition <- factor(data_correct$condition, levels = c("standing", "backward", "forward"))
data_correct$condition<-as.factor(data_correct$condition)
contrasts(data_correct$condition)

# baseline model
baseline_model<-lme(rt_offset ~ 1, random = list(~1|participant, ~1+condition|participant, ~1+Response|participant), data = data_correct, method = "ML")
summary(baseline_model)

# lme model adding Response
lme_model_Response<-update(baseline_model, .~. + Response)
anova(baseline_model, lme_model_Response) # sign. difference
summary(lme_model_Response)

# lme model adding condition
lme_model_Condition<-update(baseline_model, .~. + condition)
anova(baseline_model, lme_model_Condition) # sign. difference
summary(lme_model_Condition)

# lme model adding condition+Response
lme_model_cond_and_Resp<-update(baseline_model, .~. + Response+condition)

# lme model adding interaction condition*Response
lme_model_interaction<-update(baseline_model, .~. + Response*condition)
anova(lme_model_cond_and_Resp, lme_model_interaction) # sign. difference
summary(lme_model_interaction) 

#Check assumptions:
plot(lme_model_interaction)
qqnorm(lme_model_interaction)
# -> normal

# lme model adding order
lme_model_Order<-update(lme_model_interaction, .~. + order)
anova(lme_model_interaction, lme_model_Order) # -> n.s.

# post hoc tests 
data_backward<-subset(data_correct, data_correct$condition =="backward")
backward_post<-aggregate(rt_offset ~ Response + participant, data=data_backward, FUN=mean)
backward_post$participant<-as.factor(backward_post$participant)
pairedSamplesTTest(rt_offset ~ Response, data=backward_post, id = "participant")
# -> p = .001, Cohen's d:.59

data_forward<-subset(data_correct, data_correct$condition =="forward")
forward_post<-aggregate(rt_offset ~ Response + participant, data=data_forward, FUN=mean)
forward_post$participant<-as.factor(forward_post$participant)
pairedSamplesTTest(rt_offset ~ Response, data=forward_post, id = "participant")
# -> no differences

data_standing<-subset(data_correct, data_correct$condition =="standing")
standing_post<-aggregate(rt_offset ~ Response + participant, data=data_standing, FUN=mean)
standing_post$participant<-as.factor(standing_post$participant)
pairedSamplesTTest(rt_offset ~ Response, data=standing_post, id = "participant")
# -> no differences



